Acoustic measurement of a granular density of modes 
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In glasses and other disordered materials, measurements of the vibrational density of states reveal that an excess number of 
long-wavelength (low-frequency) modes, as compared to the Debye scaling seen in crystalline materials, is associated with a 
loss of mechanical rigidity. In this paper, we present a novel technique for measuring the density of modes (DOM) in a real 
granular material, in which we mimic thermal excitations using white noise acoustic waves. The resulting vibrations are detected 
with piezoelectric sensors embedded inside a subset of the particles, from which we are able to compute the DOM via the 
spectrum of the velocity autocorrelation function, a technique previously applied in thermal systems. The velocity distribution 
for individual particles is observed to be Gaussian, but the ensemble distribution is non-Gaussian due to varying widths of the 
individual distributions. In spite of this deviation from a true thermal system, we find that the DOM exhibits several thermal-like 
features, including Debye scaling in a compressed hexagonally ordered packing, and an increase in low-frequency modes as the 
confining pressure is decreased. In disordered packings, we find that a characteristic frequency fc increases with pressure, but 
more weakly than has been observed in simulations of frictionless packings. 
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1 Introduction 

Simulations of idealized granular materials indicate that they 
undergo a jamming transition whereby the material becomes 
rigid and able to support a finite pressure^ ^. This transition 
occurs in systems for which the average coordination num- 
ber z has increased to the critical value Zc, thereby constrain- 
ing the motion of the particles. Just above the transition, it 
has been observed in simulations of both frictionless^ ^ and 
frictionalEEl granular packings that the density of states D{cd) 
exhibits an excess number of low-frequency modes as com- 
pared to Debye scaling. Beyond some frequency (0*, the den- 
sity of states is observed to deviate from the Debye scaling 
Z)((jL)) oc 0)'^"^, where d is the dimensionality of the system. In 
frictional simulations^ ^, this crossover frequency is observed 
to scale as (0* oc (^ — ^^). The spatial eigenmodes associated 
with frequencies below 0)* are said to be soft, exhibiting long- 
wavelength rearrangements. 

While real granular materials have no thermal vibrations 
and the density of states is therefore not strictly defined, it 
has been possible to observe similar spatial modes via parti- 
cle tracking and the construction of a covariance matrix^C^. 
Such techniques are similar to those used to measure the den- 
sity of states for colloids where thermal fluctuations of 
the particles are naturally present. In experiments on an oscil- 
lated granular material, the soft modes take the form of col- 
lective rearrangements which occur as the system unjams^, 
and the number of these long- wavelength modes increases as 
the packing fraction approaches the jamming transition. A 
disadvantage of using the covariance matrix method is that it 
requires visual access to particles, and can therefore only be 
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used in highly idealized granular systems. 

In this paper, we will describe an alternative experimen- 
tal technique for measuring vibrational modes in a frictional 
granular packing. Our approach is inspired by one previously 
used to measure the density of states via the velocity auto- 
correlation function (VACF) in simulations of conventional 
solidElor li^idEI systems. In recent simulations of thermal 
soft spheres this method successfully recovered the familiar 
scaling (0* oc (0 — 0^)^^^, where is the packing fraction. The 
measurement starts from the velocity autocorrelation function 
Cy(t\ defined as 



I:,-(v,-(t + 0-v,-(t))^ 
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(1) 



where Vi(t) is the velocity of particle / as a function of time, 
(•)^ represents a temporal average, and the summation extends 
over the particles. The density of states is then given by 
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This method, which we will refer to as the VACF method, 
provides the thermal density of states, and its applicability to 
an athermal system is neither expected nor guaranteed. How- 
ever, there are several criteria which, if satisfied, would give 
us some confidence that the analogy is a reasonable one. Ide- 
ally, we want an isotropic source of vibrations which will over- 
come the dissipation and provide a steady state. Second, this 
injected energy should partition itself equally among the de- 
grees of freedom, and provide a thermal-like velocity distribu- 
tion for each particle (corresponding to the temperature of the 
system). 

In our experiments, we mimic thermal vibrations using 
acoustic excitations at the lower boundary, and test the degree 
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Fig. 1 Images of the granular packing as a function of pressure (columns) and degree of order (rows). Image (a) highlights the location of key 
parts of the apparatus. The vibrating driver is located at the bottom of the granular packing, brass weights of various sizes provide the 
confining pressure, and rough walls restrict bulk motion. Eight particles (shown by yellow circles) contain embedded piezoelectric sensors for 
detecting vibrations. Images (a,b,c) show the force chains as a function of pressure for a single disordered packing, and images (d,e,f) for a 
single ordered packing. 



to which the rest of the criteria are satisfied. Using particle- 
scale measurements from piezoelectric sensors embedded in a 
subset of particles, we apply Eq. [T] and [2] to quantify the re- 
sulting vibrational modes of the packing. While the results of 
Eq. [2] cannot properly be called a density of states (we shall 
refer to it as a density of modes), we are nonetheless able to 
investigate the empirical utility of such a quantity in describ- 
ing the state of the material near the rigidity transition. While 
the VACF method provides /)(/), it does not allow for the vi- 
sualization of the spatial modes. However, there are several 
key advantages of the method: we do not require optical ac- 
cess to the particles or need to sample all particles in the pack- 
ing. These advantages would allow our method to be easily 
adapted to a real 3D granular materials. 

2 Experimental Setup 

We perform experiments in a two dimensional granular pack- 
ing composed of discs cut from Vishay PSM-4, which is a 
photoelastic material that allows for the visualization of the 
internal force structure. The granular packing has lateral di- 
mensions of 29 X 22 cm, and is oriented vertically in order to 
both minimize friction with the walls and allow the pressure 
to be set using brass weights. In order to minimize bulk move- 
ment of the packing, the side walls of the apparatus have been 



made rough (see Fig. |T]). Further details about the apparatus 
are available in Ref. 

In this Paper, we present experiments on both ordered and 
disordered particle configurations. The ordered packings are 
composed of monodisperse discs with diameter ai = 11 mm 
arranged by hand into an ordered hexagonal pattern. The 
disordered packings are bidisperse, with an equal mixture of 
as = 9 mm and gl = 11 mm particles in order to suppress 
crystallization. The particles have density p = 1.06 g/cw?, 
thickness 6.35 mm, and static bulk modulus E = 4 kPa. We 
make measurements on a total of 60 unique particle configura- 
tions, of which 30 are disordered and 30 are ordered. For each 
of these particle configurations, we also investigate 7 different 
pressures ranging from 2.7 x 10~^ to 5.9 x 10~^ E, where the 
lowest pressure is set by the weight of the particles themselves. 

An electromagnetic driver (MB Dynamics PM50A) is at- 
tached to a driving platform of width 8.5 cm located at the 
bottom of the packing. Because a granular material is dis- 
sipative (due to both friction and viscoelasticity), continuous 
driving is necessary to maintain a steady state. We mimic ther- 
mal vibrations by subjecting the platform to oscillations which 
have a flat velocity spectrum. Observations of a piezoelectric 
sensor mounted on the driver platform indicate that the sys- 
tem (apparatus plus the granular material) has several strong 
mechanical resonances. We are able to compensate for these 
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Fig. 2 Velocity spectrum measured in a continuous sheet of PSM-4 
material, showing /"^^^ decay. 



resonances up to 3 kHz and maintain a flat driver response 
down to approximately 300 Hz. Below 300 Hz, we approx- 
imate the driver response A{f) by a third degree polynomial 
for the disordered packing and a power law for the ordered 
packing; these empirical measurements will later be used to 
correct the sensor responses. 

Vibrations are recorded from eight particles which contain a 
piezoelectric sensor at their center. The particles are arranged 
in two rows, the first one 4.5 cm from the bottom wall and the 
second at 10 cm (see Fig.[T]). Piezoelectric materials are sensi- 
tive to a single direction of stress, and therefore only measure 
the component of acceleration along that axis. During setup, 
individual sensors are placed in a variety of orientations. This 
ensures that the ensemble of sensors collects data from many 
different directions, both multiply- scattered waves and some 
which arrive directly from the driver. 

Piezoelectric sensors produce a voltage V{t) proportional to 
the stress on the sensor, and therefore also proportional to its 
instantaneous acceleration. Each sensor is sampled for a du- 
ration of 2 sec at a sampling rate of 100 kHz. By integrating 
v'{t) = JV{t)dt, we obtain a velocity measurement v\ which 
must still be corrected by deconvolving both the viscoelas- 
tic response of the bulk material and the low-frequency driver 
response A{f). As the measured voltages are difficult to cali- 
brate, all velocities will be reported in arbitrary, but consistent, 
units. 

In order to factor out the viscoelasticity of the photoelastic 
material, we measure the response of a continuous sheet of 
PSM-4 embedded with identical piezoelectric sensors, subject 
to the same driving input as the granular experiments. The 
result of this test is shown in Fig. [2] We observe that the ve- 
locity spectrum exhibits a f~^^^ decay. Strong, frequency de- 
pendent, damping is expected for viscoelastic materials'^. 



In Fourier space (denoted by ^), we combine these mea- 
surements to obtain v = v'/^/^/A. We additionally apply a 
band-pass filter with a low-frequency cutoff of / = 75Hz and 
a high-frequency cutoff of / = 3 kHz, as these regimes are 
dominated by electrical noise and the mechanical response of 
the apparatus, respectively. Finally, we perform an inverse 
transform and use the resulting v(^) for the remainder of the 
analysis. 

By using birefringent particles, we can measure the parti- 
cle positions via a Hough transform and the vector contact 
forces via a nonlinear fit to the photo elastic fringes visible 
through crossed circular polarizersE^ESl. We make these mea- 
surements on particles located in the same region as the sensor 
particles. We exclude rattlers (particles with less than two de- 
tectable neighbors) from the analysis. 

In order to make comparisons with computer simulations, 
it is helpful to non-dimensionalize our measurements. Fol- 
lowing conventional, we take our frequency unit to be /o = 

7y-\/ where m is the mass of the particles, a is the mean 
particle diameter, and Vq comes from the harmonic particle 
contact law used in the simulation. To estimate Vq for our 
particles, we compare the force law that corresponds to the 
interaction potential (F = ^5, where 5 is the Hertzian over- 

TlEL ^ 



lap) to the force law for ideal Hertzian discs (F 
where L is the particle thickness and V is the Poisson ratio). 
Matching the coefficients, we find that /o(/) = 



8(1- 



LEjf) 



647rm(l-v2) • 

This treatment neglects the observation that our particles de- 
viate from the linear Hertzian law, and instead have a force 
law F oc 5^/^. Due to the viscoelasticity of our particles, the 
modulus is frequency-dependent, with modulus E{f) oc /V2 
Therefore, /o is itself frequency-dependent, and we report fre- 
quencies as the ratio in order to make comparisons with 
simulations. The range of accessible non-dimensional fre- 
quencies corresponds to 0.06 to 0.8. 

3 Results 

3.1 Thermal analogy 

First, we characterize the extent to which the analogy to ther- 
mal systems is reasonable. We anisotropically inject energy 
from a single boundary of the granular packing for 2.8 sec, of 
which the middle 2 sec is analyzed. The typical time for V{t) 
measurements to either reach steady state or to decay when 
driving ceases is on the order of 1 ms. 

In a thermal solid, particles should have a Gaussian velocity 
distribution, where the standard deviation of the distribution is 
set by the temperature of the material. For each sensored parti- 
cle, in each experiment, we measure the velocity distributions 
<^i{v) (see Fig.[3|, and find that these individual distributions 
are well-described by Gaussian distributions; the average kur- 
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Fig. 3 (a) Two velocity distributions ^i{v), for two sensored particles in the first row of the same disordered packing, with pressure 

P = 2.2 X 10~^ E. Ensemble velocity distributions ^(v) of all sensored particles, rows, and experiments in (b) disordered packings and (c) 

ordered packings. Gaussian distributions with the same standard deviation as the data are shown as solid lines on all plots. 
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Fig. 4 Example velocity autocorrelation function (from Eq.[lJ, in a 
disordered packing with P = 5.9x10"^ E. Correlations decay within 
approximately 20 ms for all experimental runs. 



tosis of the distributions is 3.0 ±0.07. However, different par- 
ticles, even from the same row of the same experiment, have 
different widths, depending on the local environment. From 
the measured vector contact forces on each sensored particle, 
we observe that the width of a particle's velocity distribution 
increases with the local pressure from its neighbors. The trend 
is consistent with what was previously observed for mea- 
surements of sound wave amplitudes, and can be attributed 
to an increased transmission of energy through contacts with 
larger contact area. 

As a result of these spatial heterogeneities, the velocity dis- 
tribution ^(v) measured over all 30 experiments at the same 
P is not a Gaussian. Instead, as shown in Fig.|3})c, the tails are 
significantly broader; this holds for both the disordered and 
ordered packings. In addition, the ensemble distributions have 
a width which increases with increasing pressure. This means 



that the ensemble does not satisfy equipartition, and we must 
therefore proceed with caution. 

3.2 Determining D{f) 

In order to utilize the VACF method, we would ideally use 
all particles in a packing when calculating the velocity auto- 
correlation Cy{t) from Eq. [T] Instead, we average over eight 
particles within 30 experiments at the same pressure and same 
order/disorder, and obtain a long time-average for each. An 
example Cy{t) is shown in Fig.jij note that the decay time is 
much shorter (20 ms) than the 2 sec measurement duration, 
providing many decorrelated measurements during each ex- 
periment. Using the real part of the spectrum of Cv{t), we cal- 
culate the density of modes D{f) via Eq.[2] There are strong 
peaks at multiples of 60 Hz due to electrical noise which we 
filter out as a last step after calculating D{f). 

We calculate D{f) for each of seven confining pressures P, 
applied to 30 ordered and 30 disordered packings. The results 
are shown in Fig. [s^b. We observe Debye scaling (/^~^ for 
dimension d = 2) for large P in the ordered packings, which 
is expected since these packings are the most uniform in the 
coordination number z and and contact forces. Debye scaling 
is not observed in the disordered experiments. 

In order to compare the ordered and disordered packings, it 
is helpful to consider the ratio D^is/Dord- In Fig.js]:, we plot 
this ratio and thereby eliminate some of the remaining elec- 
tronic/apparatus resonances: note that the sharp decrease near 
1500 Hz, as well as number of smaller features, are removed. 
We observe from this ratio that at all values of P, disordered 
packings have more low-frequency modes than do the ordered 
packings. 

As a function of decreasing P, both ordered and disordered 
packings exhibit a growing number of low-frequency modes. 
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Fig. 5 Measured density of modes D{f) as a function of confining pressure P for (a) ordered and (b) disordered packings. Line color 
corresponds to P as shown in the legend. The frequencies are reported in both experimental frequency units (bottom, Hz) and simulation units 
(top, ///o(/))- In (a), the solid black line provides a comparison to Debye scaling for 2D systems. In (b), the solid line shows the threshold 
specifying the cutoff frequency fc. (c) Ratio of the ordered and disordered D{f). (d) Scaling of fc as a function of P (bottom axis) and z — 3 
(top axis). 



In the case of the ordered packings, it is important to note that 
at low P, disorder is still present in both the contact network 
(z < 6 for many of the particles) and the force chain network 
(see Fig. [T]). For the disordered packings, we quantify this 
decrease in the number of low-frequency modes by calculating 
a critical frequency fc at which D{f) falls below a threshold 
(see Fig. [5])). 

Fig.[5]l shows how fc depends on the confining pressure P\ 
we fit a power law fc ^ p0 ii±0 05^ where the mean and un- 
certainty is determined from a bootstrapping method that re- 
peats the analysis for subsets of the data, randomly resampled 
with replacement. For a frictionless packing with Hertzian ex- 
ponent a = 5/4, /c would be expected to scale as fc ^ P^, 
with an exponent j3 =2/5; our observations show a consid- 
erably flatter trend. This can be understood in light of simu- 
lations of frictional particles^, where power-law scaling was 
only observed in the large-friction limit. For experimentally- 
reasonable friction coefficients ({l = 0.2 to 0.8), no simple 
(fc,P) scaling was observed. However, they did observe a 
friction-independent scaling relationship cOc (z — 3) (where 
Zc = 3 is the isostatic value for a frictional packing); Henkes 
et al.^ observed a similar scaling relation. In our experi- 
mental results, good estimates of z — 3 are only available for 
P > 10~^E; this data is plotted in Fig.[5]l for comparison and 
shows the expected steeper trend. The range of measured val- 
ues is unfortunately too small to draw a more quantitative con- 
clusion. 

4 Discussion 

In calculating a density of modes D{f) for our athermal, fric- 
tional, viscoelastic, granular packings, we mimic thermal- 
ization using acoustic excitation. We find that while some 
thermal-like conditions are satisfied - we can achieve an 
isotropic, steady state through continuous driving and mea- 
suring vibrations resulting from multiply- scattered waves - 



the velocity distributions are fundamentally non-thermal-like. 
While individual particles have Gaussian velocity distribu- 
tions, each one is sufficiently different that the ensemble dis- 
tribution does not. Notably, experiments by Brito et al.^ 
also found that displacement distributions were individually 
Gaussian, but with widths which varied from particle to par- 
ticle. Through the use of photoelastic particles, we explain 
these variations as arising from spatial heterogeneities in the 
force chains. In spite of these non-thermal features, the VACF 
method provides a means to calculate D{f), and the observa- 
tion of Debye scaling for the most uniform (well-compressed, 
ordered) packings suggest that the analogy is justified. 

The VACF method differs significantly from what has been 
utilized in prior experimentslSnM which particle trajecto- 
ries are used to construct a co variance matrix The covari- 
ance matrix method works well in systems where there is opti- 
cal access to accurate position information for a large number 
of particles, and has the benefit of providing a means to visu- 
alize the spatial modes at each frequency. However, real gran- 
ular packings do not typically provide the necessary optical 
access for such methods to be feasible. In addition, the quality 
of displacement measurements depends strongly on the frame 
rate of the camera, which is often at odds with the high spatial 
resolution needed for accurate particle-tracking. In contrast, 
our application of the VACF method does not require optical 
access to the packing, and instead relies on only a small num- 
ber of piezoelectric sensors which readily provide high tem- 
poral resolution. While we have used photoelastic particles in 
order to visualize the internal stresses, the VACF method does 
not require such a packing, and would be applicable for use in 
ordinary, non-circular, three-dimensional granular packings. 

A useful application of a granular D{f) would be to mea- 
sure how far a packing is from P = 0, ^c^ or Zc- While we 
observe only a flat trend {fc ^ P^-^^), there is evidence for a 
stronger dependence on the coordination number. This is very 
encouraging, since the presence/absence of contacts has been 
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extremely difficult to measure accurately in experiments. Fur- 
ther investigations should explore this scaling relation more 
quantitatively. 

5 Conclusion 

Our experiments use a novel technique to apply thermal meth- 
ods to a real, athermal, granular packing through the use of 
acoustic excitation. We have probed the mode structure of 
a granular packing by acoustically mimicking thermal vibra- 
tions and calculated a density of modes D{f) using the spec- 
trum of the velocity autocorrelation function. Our granular 
packings deviate in many important ways from the ideal ones 
used in simulations since our particles are soft, frictional, and 
viscoelastic; furthermore, there is no equipartition of vibra- 
tion amplitudes among the particles. Nonetheless, we can re- 
cover many important features such as Debye scaling at high 
pressure in the ordered packings and an excess number of low- 
frequency modes compared to Debye scaling in the disordered 
packings. These low-frequency modes can be used to differ- 
entiate the pressure state of the packing as it is seen that a 
crossover frequency fc varies with both pressure and coordi- 
nation number. We have thereby been able to use a small num- 
ber of particle scale measurements to see system-level proper- 
ties, a technique which could be extended to granular materials 
where there is no optical access. 
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